Added xt::convolve by spectre-ns · Pull Request #2507 · xtensor-stack/xtensor

spectre-ns · 2022-04-03T02:32:36Z

Checklist

The title and commit message(s) are descriptive.
Small commits made to fix your PR have been squashed to avoid history pollution.
Tests have been added for new features or bug fixes.
API of new functions and classes are documented.

Description

Added a 1D convolution function with a similar API as numpy and provide a couple small test to verify it's functionality.
Fixes #1928

spectre-ns · 2022-04-03T02:34:15Z

I noted where the code could be improved with new constexpr functionality released in C++17 I think this library strives for C++14 or newer compatibility correct?

JohanMabille

Thanks for tackling this!

JohanMabille · 2022-04-05T07:29:31Z

+    *   should v be longer than a.
+    */
+    template <class E1, class E2, class E3>
+    inline auto convolve(E1&& a, E2&& v, E3&& mode)


E3 being the mode, it can be passed by value.

JohanMabille · 2022-04-05T07:34:33Z

+            auto temp = a;
+            a = v;
+            v = temp;
+        }


A mean to avoid copy is to have an additional function performing the implementation that is called by this one:

namespace detail { template <class E1, class E2, class E3> inline auto convolve_impl(E1&& e1, E2&& e2, E3 mode) { // Implementation with no checks } } template <class E1, class E2, class E3> inline auto convolve(E1&& a, E2&& b, E3 mode) { // Dimensions and shape checks // ... bool cond = a.shape(0) < v.shape(0); return cond ? detail::convole_impl(std::forward<E2>(v), std::forward<E1>(a), mode) : detail::convole_impl(std::forward<E1>(a), std::forward<E2>(v), mode) }

JohanMabille · 2022-04-05T07:35:36Z

+    namespace conv_mode
+    {
+        struct valid{};
+        struct full {};


Is that intentional to not implement the same mode?

Yeah I figured I could push that implementation off to another pull request because it can be implemented by trimming the full output.

JohanMabille · 2022-04-05T07:39:29Z

+            v = temp;
+        }
+
+        //this could benefit from C++17 if constexpr


This can be achieved with tag dispatching (although I'm pretty sure the compiler can optimize the current implementation)

template <class E1, class E2, class E3> inline auto convolve_impl(E1&& e1, E2&& e2, E3 mode) { XTENSOR_THROW(...); } template <class E1, class E2> inline auto convolve_impl(E1&& e1, E2&& e2, conv_mode::full) { // ... } template <class E1, class E2> inline auto convolve_impl(E1&& e1, E2&& e2, conv_mode::valid) { //.... }

Done! Didn't implement the generic case as it would be better to have compile time error over runtime.

JohanMabille · 2022-04-05T07:41:55Z

+            size_t const nv = v.shape(0);
+            size_t const n = na + nv - 1;
+            xt::xarray<value_type> out = xt::zeros<value_type>({ n });
+            for (size_t i = 0; i < n; ++i) {


nitpicking: an opening curly brace goes on a new line

tdegeus

Awesome!

tdegeus · 2022-04-06T08:54:37Z

+    /*
+    * convolution mode placeholders
+    */
+    namespace conv_mode


I would prefer the name convolve_mode I think.

tdegeus · 2022-04-06T08:57:07Z

+    * @param a 1D expression
+    * @param v 1D convolutional kernel


Very minor: This is a bit different than NumPy's documentation. In principle that is ok, it's just that one could argue that either of the arguments is the kernel of the other. So naming one 'kernel' might confuse?

This is a non-issue no because I have removed the copy will make it symmetrical

tdegeus · 2022-04-06T08:58:46Z

+    *
+    * @param a 1D expression
+    * @param v 1D convolutional kernel
+    * @param mode placeholder from conv_mode namespace


Suggested change

* @param mode placeholder from conv_mode namespace

* @param mode placeholder Select algorithm from #conv_mode

tdegeus · 2022-04-06T08:59:29Z

    }
+
+    /*
+    * convolution mode placeholders


A brief description is on order I think. Also you could point out that these modes are the same as for NumPy.

tdegeus · 2022-04-06T09:00:47Z

+    {
+        using value_type = typename std::decay<E1>::type::value_type;
+
+        if (a.shape().size() != 1 || v.shape().size() != 1)


Suggested change

if (a.shape().size() != 1 || v.shape().size() != 1)

if (a.dimension() != 1 || v.dimension() != 1)

tdegeus · 2022-04-06T09:06:02Z

+                size_t const jmn = (i >= nv - 1) ? i - (nv - 1) : 0;
+                size_t const jmx = (i < na - 1) ? i : na - 1;
+                for (size_t j = jmn; j <= jmx; ++j) {
+                    out[i] += (a(j) * v(i - j));


Suggested change

out[i] += (a(j) * v(i - j));

out(i) += a(j) * v(i - j);

(I don't think the braces are needed there

tdegeus · 2022-04-06T09:06:42Z

+            size_t const na = a.shape(0);
+            size_t const nv = v.shape(0);
+            int const n = std::max(na, nv) - std::min(na, nv) + 1;
+            xt::xarray<value_type> out = xt::zeros<value_type>({ n });


tdegeus · 2022-04-06T09:07:18Z

+            xt::xarray<value_type> out = xt::zeros<value_type>({ n });
+            for (size_t i = 0; i < n; ++i) {
+                for (int j(v.shape(0) - 1), k(i); j >= 0; --j) {
+                    out[i] += a(j) * v(k);


Suggested change

out[i] += a(j) * v(k);

out(i) += a(j) * v(k);

tdegeus · 2022-04-06T09:11:51Z

+            int const n = std::max(na, nv) - std::min(na, nv) + 1;
+            xt::xarray<value_type> out = xt::zeros<value_type>({ n });
+            for (size_t i = 0; i < n; ++i) {
+                for (int j(v.shape(0) - 1), k(i); j >= 0; --j) {


Suggested change

for (int j(v.shape(0) - 1), k(i); j >= 0; --j) {

for (int j = v.shape(0) - 1, k = i; j >= 0; --j) {

no?
Also I think that the int can be avoided, e.g. https://stackoverflow.com/a/23244694/2646505

But most importantly, why is it even relevant to loop backward?

tdegeus · 2022-04-06T09:12:44Z

+
+    TEST(xmath, convolve_full)
+    {
+        xt::xarray<double> x = { 1.0, 1.0, 1.0 };


It be'd good to add another test which does not only involve unity inputs

tdegeus

Thanks! Did not have time to check the conditions in detail

spectre-ns · 2022-04-12T22:28:55Z

@tdegeus I accepted all your suggestions! Thanks!

tdegeus

Nitpicking. Other then that things look good to me, thanks! A side note: I'm surprised that NumPy does not support a periodic version. True, one could (and probably should) use Fourier for that. But still

spectre-ns · 2022-04-13T10:36:27Z

Nitpicking. Other then that things look good to me, thanks! A side note: I'm surprised that NumPy does not support a periodic version. True, one could (and probably should) use Fourier for that. But still

I think a periodic convolution in the time domain would simply be a wrapped convolution? But yeah Fourier does the wrapping implicitly so best to do that. Scipy probably offers that functionality.

tdegeus · 2022-04-13T12:02:44Z

Final question. Could you squash your commits? We prefer here to not have GitHub do that, such that we have an easier way to generate a changeling (if I remember well, I think that nowadays GitHub could do it, still in the mean-time...). Thanks!

spectre-ns · 2022-04-13T12:06:57Z

Final question. Could you squash your commits? We prefer here to not have GitHub do that, such that we have an easier way to generate a changeling (if I remember well, I think that nowadays GitHub could do it, still in the mean-time...). Thanks!

I definitely will do that. When I submitted there were only two commits. I appreciate the thorough review!

spectre-ns · 2022-04-13T12:45:40Z

This is embarrassing... I messed hopefully I didn't break the pr...

tdegeus · 2022-04-13T12:56:48Z

Still looks ok, the merge commit was not necessary I guess.
It's just empty so squashing it again should get rid of it.

spectre-ns · 2022-04-13T20:01:43Z

Still looks ok, the merge commit was not necessary I guess. It's just empty so squashing it again should get rid of it.

Fixed it

JohanMabille · 2022-04-15T07:14:36Z

Awesome! Thanks a lot!

JohanMabille requested changes Apr 5, 2022

View reviewed changes

tdegeus requested changes Apr 6, 2022

View reviewed changes

spectre-ns requested review from JohanMabille and tdegeus April 10, 2022 17:49

tdegeus requested changes Apr 11, 2022

View reviewed changes

spectre-ns requested a review from tdegeus April 12, 2022 22:29

tdegeus reviewed Apr 13, 2022

View reviewed changes

Comment thread include/xtensor/xmath.hpp Outdated

Comment thread include/xtensor/xmath.hpp Outdated

tdegeus approved these changes Apr 13, 2022

View reviewed changes

spectre-ns closed this Apr 13, 2022

spectre-ns reopened this Apr 13, 2022

spectre-ns closed this Apr 13, 2022

spectre-ns reopened this Apr 13, 2022

spectre-ns closed this Apr 13, 2022

Added xt::convolve

935c649

spectre-ns reopened this Apr 13, 2022

JohanMabille approved these changes Apr 15, 2022

View reviewed changes

JohanMabille merged commit fd7a610 into xtensor-stack:master Apr 15, 2022

	* @param mode placeholder from conv_mode namespace
	* @param mode placeholder Select algorithm from #conv_mode

	if (a.shape().size() != 1 \|\| v.shape().size() != 1)
	if (a.dimension() != 1 \|\| v.dimension() != 1)

	for (int j(v.shape(0) - 1), k(i); j >= 0; --j) {
	for (int j = v.shape(0) - 1, k = i; j >= 0; --j) {

Conversation

spectre-ns commented Apr 3, 2022 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Checklist

Description

Uh oh!

spectre-ns commented Apr 3, 2022

Uh oh!

JohanMabille left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

JohanMabille Apr 5, 2022 • edited by tdegeus Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

tdegeus left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

tdegeus left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

spectre-ns commented Apr 12, 2022

Uh oh!

tdegeus left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Uh oh!

spectre-ns commented Apr 13, 2022

Uh oh!

tdegeus commented Apr 13, 2022

Uh oh!

spectre-ns commented Apr 13, 2022

Uh oh!

spectre-ns commented Apr 3, 2022 •

edited

Loading

JohanMabille Apr 5, 2022 •

edited by tdegeus

Loading

spectre-ns commented Apr 13, 2022 •

edited

Loading

tdegeus commented Apr 13, 2022 •

edited

Loading